AO- A 1 36  387 
UNCLASSIFIED 


DYNAMIC  ANALYSIS  OF  STRUCTURES  USING  LANCZOS 
COORDINATES(U)  CALIFORNIA  UNIV  BERKELEY  CENTER  FOR  PURE 
AND  APPLIED  MATHEMATICS  B  NOUR-OMID  ET  AL .  NOV  83 
PAM  186  N000I4-76-C-0013  F/G  12/1 


MICROCOPY  RESOLUTION  TEST  CHART 

NATIONAL  BUREAU  OF  STANDARDS- 1963- A 


DYNAMIC  ANALYSIS  OF  STRUCTURES  USING 


LANCZOS  COORDINATES 
by 

Bahrain  Nour-Omid  * 
and 

Ray  W.  Clough 1 

SUMMARY 

A  procedure  for  deriving  the  Lancsou  Coordinates  is  explained  and  their  use  in  structural 
dynamics  analysis  as  an  alternative  to  modal  coordinates  is  discussed.  The  coordinates  are 
obtained  by  an  inverse  iteration  procedure  in  which  orthogonality  is  imposed  between  the  vectors 
resulting  from  successive  iteration  cycles.  Using  these  Lancxos  coordinates  the  equations  of 
motion  are  transformed  to  tridiagonal  form,  which  provides  for  very  efficient  time-stepping  solu¬ 
tion.  The  effectiveness  of  the  method  is  demonstrated  by  a  numerical  example. 
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Introdoctlon 

Analysis  of  the  response  of  structures  subjected  to  dynnmie  loads  typically  involves  first 
establishing  a  discretized  model  of  the  structural  system,  using  an  appropriate  discretization  pro¬ 
cedure  such  as  finite  elements  or  finite  differences.  Then  the  equations  of  motion  are  expressed  as 
follows: 

Mi!  +  Cu  +  Ku  -  f  (1) 

in  which  M.  C,  K,  are  the  n  X  n  mass,  damping  and  stiffness  matrices,  respectively;  f  is  the  time 
dependent  external  load;  and  u  is  the  displacement  vector  describing  the  response  of  the  struc¬ 
ture.  The  number  of  displacement  coordinates,  n ,  employed  in  the  discretization  depends  on  the 
complexity  of  the  structure  and  on  the  amount  of  detail  desired  in  the  description  of  the  dynamic 
stress  response. 

Because  the  discretized  model  of  a  complicated  structural  system  may  have  many  hundreds 
or  even  thousands  of  degrees  of  freedom,  it  is  customary  to  reduce  the  equations  of  motion  to  a 
much  smaller  number  before  the  dynamic  response  is  calculated.  In  the  analysis  of  linear  struc¬ 
tures,  it  has  become  essentially  standard  practice  to  express  the  response  in  terms  of  the 
undamped  free  vibration  mode  shapes,  using  only  enough  of  the  lower  modes  to  express  the 
behavior  adequately.  The  main  analytical  problem  then  becomes  the  evaluation  of  the  mode 
shapes,  and  the  problem  of  reducing  the  number  of  degree  of  freedom  is  transferred  to  this  phase 
of  the  analysis. 

The  Rayleigh-Ritz  method  |1)  has  been  used  widely  to  reduce  the  dimension  of  the  equa¬ 
tions  of  motion.  This  may  reduce  the  dimensions  of  the  equations  of  motion.  This  may  be  looked 
upon  as  nothing  more  than  a  coordinate  transformation,  the  essential  step  being  the  choice  of  a 
set  of  vectors  Y„  “  IXi.ys.ya.  *'*•?«]  f<*  describing  the  dynamic  response.  Thus  an 
approximation  to  the  dynamic  response,  ua,  is  expressed  as  a  linear  combination  of  the  chosen 
displacement  vectors 

(2) 


i 
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ii  which  the  amplitudes  of  these  vectors,  x(t),  are  the  generalized  coordinates  representing  the 
approximate  solution.  Applying  this  coordinate  transformation  and  its  time  derivatives  to  Eq.  1, 
leads  to 

Sf„u  +  C„u  +  Ku u  -  f.  (3) 

in  which  the  generalized  coordinate  mass,  damping,  and  stiffness  matrices  and  the  external  force 
vector  are  defined  as  follows: 

-  yXmy. 

-YlCY, 

k«  -  y;ky„  (4) 

K  -Y't 

The  bask  problem  in  the  use  of  the  Rayleigh-Ritz  method  is  that  the  quality  of  the  approx¬ 
imate  solution,  u„ ,  depends  entirely  on  the  adequacy  of  the  assumed  set  of  displacement  vectors, 
Y„.  One  approach  to  dealing  with  this  problem  is  to  apply  iterative  “improvement”  of  the  vec¬ 
tors.  In  this  method,  known  as  “subspace”  or  “simultaneous”  iteration,  the  trial  vectors  are  all 
subjected  to  inverse  iteration  combined  with  some  technique  (such  as  Gram-Schmidt  orthogonali- 
zation)  that  forces  convergence  to  independent  shapes.  Actually  the  convergence  is  to  the  lowest 
undamped  vibration  mode  shapes,  §m  and  then  the  coordinate  transformation 

(5) 

leads  to  an  uncoupled  set  of  modal  equations.  That  is,  the  equations  equivalent  to  Eq.  3  are 
independent  because  the  mode  shape  orthogonality  properties  applied  to  Eq.  4  leads  to  diagonal 
generalized  mass,  damping  and  stiffness  matrices  (assuming  C  is  a  proportional  damping  matrix). 

It  frequently  is  assumed  that  this  modal  coordinate  transformation  provides  the  most 
efficient  method  of  dynamic  response  analysis  because  the  independent  modal  equations  can  be 
solved  separately  and  the  total  response  obtained  by  superposition.  However,  there  is  no  reason 
to  believe  that  the  modal  coordinates  will  always  give  the  best  results  with  the  fewest  degrees  of 
freedom.  In  fact,  alternatives  or  supplements  to  mode  superposition,  such  as  the  mode  accelera¬ 
tion  methods  [2j  and  static  corrections  (3],  have  been  developed  to  obtain  satisfactory  results  with 
fewer  coordinates.  It  m  the  purpose  of  this  paper  to  describe  how  Lanczos  coordinates  can  lead  to 
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a  still  more  efficient  procedure  for  dynamic  response  analysis. 

Lane  sos  Algorithm 

The  Lanczos  method  [4|,  first  introduced  in  1950,  is  an  efficient  algorithm  for  extracting 
some  frequencies  and  mode  shapes  of  an  eigensystem.  It  can  be  thought  of  as  a  means  of  con¬ 
structing  an  orthogonal  set  of  vectors,  known  as  Lanczos  vectors,  for  use  in  the  Rayleigh-Ritz 
|l,7|  procedure.  Recently  Wilson,  et  al  [5]  independently  developed  the  idea  of  generating  a  set  of 
orthogonal  vectors  for  use  in  a  Ritz  type  analysis  as  an  alternative  to  the  mode  superposition 
method.  Their  method  is  identical  to  the  Lanczos  algorithm  with  full  reorthogonalization  [6,7]. 
The  results  they  report  show  that  these  “Ritz  vectors"  give  better  results  than  the  same  number 
of  modal  coordinates. 

It  will  be  shown  here  that  the  Lanczos  vectors  can  be  used  to  formulate  a  very  efficient 
means  of  dynamic  response  analysis  without  solving  the  vibration  eigenproblem.  These  vectors  do 
not  have  the  full  uncoupling  property  of  the  mode  shapes,  but  they  are  much  less  expensive  to 
generate.  Moreover,  they  are  derived  from  the  applied  load  vector  and  include  the  static  dis¬ 
placed  shape  as  the  first  vector;  thus  no  “static  correction"  is  required  no  matter  what  the  shape 
of  the  loading  is  nor  how  few  vectors  are  employed  in  the  analysis.  Moreover,  an  indication  of  the 
number  of  vectors  required  to  obtain  the  desired  degree  of  accuracy  may  be  determined  during 
the  derivation  of  the  Lanczos  vectors. 

The  Lanczos  algorithm  is  closely  related  to  the  inverse  iteration  and  power  methods  for  cal¬ 
culating  a  single  vibration  mode  shape  and  frequency  of  an  eigen  problem.  Given  a  pair  of 
matrices  K  and  M,  and  a  starting  vector  r  these  methods  generates  a  sequence  of  vectors, 

[  r ,  K_I  M  r ,  (K*‘  M j*r . (K_I  M)'  r  [,  during  j  iterations.  These  vectors  are  referred  to  as  the 

Krylov  sequence;  the  sequence  converges  to  the  eigenvector  corresponding  to  the  smallest  eigen¬ 
value,  in  magnitude,  of  the  eigenproblem  (K  -  XM)s  —  0. 

The  basic  difference  between  the  Lanczos  method  and  the  other  two  methods  is  that  the 
information  contained  in  each  successive  vector  of  the  Krylov  sequence  is  employed  in  the 
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dynamic  response  analysis,  instead  of  using  only  the  Baal  converged  vector.  To  be  more  specific, 
the  Laacxos  algorithm  involves  supplementing  the  Krylov  sequence  with  a  Gram-Schmidt  orthog- 
onalisatioa  process  at  each  step;  the  result  is  a  set  of  M-orthonormal  vectors  that  is  used  to 
reduce  the  dimension  of  the  dynamic  equation  set.  In  the  form  of  the  Lancxos  algorithm 
described  here,  orthogonalisation  is  applied  only  with  respect  to  two  preceding  vectors,  leading  to 
a  tridiagonal  form  of  the  dynamic  equations  that  can  be  used  to  great  advantage  either  directly  in 
time  step  integration  or  in  solution  of  the  free  vibration  eigenproblem. 

To  derive  the  Lancxos  algorithm  it  will  be  assumed  for  the  moment  that  the  first  j  Lancxos 
vectors,  (<li ,  <fc ,  •  ■  •  ,  q, )  have  been  found,  and  the  analysis  of  the  j  +  1  vector  will  be  performed. 
The  resulting  vectors  all  satisfy  the  condition  q,rMqj  »  St}  where  St)  is  the  Kroneker  delta;  that 
is  the  vectors  are  orthonormal  with  respect  to  the  mass  matrix.  To  calculate  q;+l,  a  preliminary 
vector  F,  is  first  calculated  from  the  previous  vector,  q, ,  as  in  the  Krylov  sequence: 


F,  -  K  ‘Mq;  (6) 

Now  in  general  it  may  be  assumed  that  this  preliminary  vector  contains  components  from  each  of 
the  preceding  vectors,  thus 

F,  -  r;  +  a,q;  +  0,q,_,  +  7,q,-s  +  *  '  *  (7) 

where  r,  is  the  “pure”  vector  orthogonal  to  all  previous  Lanczoe  vectors,  and  a;,  f)Jt  y;, 
are  the  amplitudes  of  the  previous  vectors  contained  in  F;.  These  amplitude  coefficients  are 
evaluated  from  the  ortbonormality  of  the  Lanczoe  vectors.  Thus  if  both  sides  of  Eq.  7  are  multi* 
plied  by  q/M,  the  result  is 

q/Mr,  —  q/M t,  +  a,q*Mq,  +  £;q/Mq,-t  +  7;q/fcdq,_2  +  •  •  ■  (8) 

Here  the  first  term  on  the  right  hand  side  vanishes  by  definition  and  all  terms  beyond  the  second 
vanish  similarly  due  to  orthogonality.  The  normalizing  definition  applied  to  the  second  term  then 
reduces  Eq.  8  to  an  expression  for  the  amplitude  of  q;  along  F;: 

The  amplitude  of  q,_,  contained  in  F;  may  be  found  similarly  by  multiplying  Eq.  7  by 
q/iM.  In  this  case  all  terms  except  the  third  vanish  by  orthogonality  and  the  coefficient  of  0,  is 
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unity,  so 

$}  -  q;-iMF, 

But  substituting  Eq.  6  this  gives  0,  »  and  applying  the  transpose  of  Eq.  6  to  the 

q£i  vector  gives 

0,  -  F/-iMq,  (10) 

Finally,  expanding  F;_i  in  terms  of  its  pure  vector,  r;_i,  and  the  preceding  Lanczos  vectors,  as  in 
Eq.  7,  the  transpose  of  Eq.  10  becomes 

0,  —  q/Mr;-i  +  o;_,q/Mq,-i  +  0}.1q)TMq}.i  +  7_,-iq/fcfq7-3  +  '  '  '  (11) 

It  is  evident  that  all  terms  except  the  first  vanish  on  the  right  hand  side.  Now  q,  is  the  vector 
obtained  by  normalizing  r,_t  with  0,  as  the  normalizing  factor,  i.e. 


_1_ 

0, 

then  the  value  of  0,  is  given  by  Eq.  11  :  0,  —  or 

0) 


(12) 


0*  -  r^Mr,.,  (13) 

Continuing  in  the  same  way,  the  amplitude  of  q,_j  contained  in  F;  is  found  to  be 

7,  —  q^-sMF,  (14) 

Following  the  procedure  used  to  derive  Eq.  11,  this  leads  to 


7,  —  q/fc**,-8  +  a^flAq^  +  0,-rtfMq,^  +  7,-sq/k*q,-4  +  •  •  •  (15) 

But  using  the  normalizing  relationship  equivalent  to  Eq.  12,  r;.2  =  0,.\q,.\,  hence  when  this  is 
substituted  into  Eq.  15  it  is  clear  that  all  terms  on  the  right  hand  side  vanish,  with  the  result  that 
1j  *  0.  A  corresponding  procedure  could  be  used  to  demonstrate  that  all  further  terms  in  the 
expansion  for  F;  (Eq.  7)  vanish;  in  other  words,  the  orthogonalization  procedure  used  in  generat¬ 
ing  each  Lanczos  vector  need  be  applied  only  to  the  previous  two  vectors. 

In  summary,  the  Lanczos  algorithm  to  derive  the  vector  q,+l,  may  be  expressed  by  the  fol¬ 
lowing  sequence  of  equations: 
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T,  -  K*‘Mq, 

(») 

r;  —  F,  -  <*;q;  - 

(b) 

where  :  a}  «■  qfh tF} 

(*) 

0,  —  (r£iMr,-i) 

(<0 

1 "  XX’’ 

(e) 

where  :  0,  »  (rjrMr] ) 

(f) 

Starting  Vector 

In  general,  the  starting  vector  q,  for  the  Lanezos  algorithm  may  be  chosen  arbitrarily. 
Then  setting  q,sO,  the  procedure  defined  by  Eqs.  16  will  produce  the  second  Lanezos  vector 
and  the  process  may  be  continued  to  produce  as  many  vectors  as  are  desired.  It  will  be  noted 
that  the  cost  of  evaluating  each  successive  vector  is  constant  after  tU  has  been  obtained  because  it 
is  necessary  to  impose  orthogonality  only  with  respect  to  the  two  preceding  vectors. 

The  Lanezos  algorithm  is  particularly  advantageous  in  dynamic  analysis  when  the  applied 
load  is  of  the  form 

f  —  bc(<)  (17) 

that  is,  when  the  distribution  of  the  load,  b  remains  constant  and  only  its  amplitude  varies  with 
time.  The  starting  vector  in  this  case  is  taken  in  the  direction  of  the  static  displacement,  thus 

r0  -  /91ql  -  K‘lb  (18) 

This  choice  gives  the  Lanezos  vectors  the  important  advantage  noted  before,  that  they  automati¬ 
cally  include  the  static  displacement  and  avoid  any  possible  need  for  a  static  correction.  It  is  of 
interest  to  note  that  Wilson,  et  al.  (5j  adopted  this  starting  vector  for  this  reason.  Further,  it  is 
evident  from  Eqi.  16  that  if  b  is  orthogonal  to  any  of  the  modes  then  all  the  Lanezos  vectors  will 
also  be  orthogonal  to  this  mode. 

An  important  type  of  dynamic  loading  of  the  form  of  Eq.  17  is  earthquake  excitation,  where 
the  effective  loading  results  from  accelerations  applied  at  the  supports.  In  this  case,  b  *  Ms 
where  a  is  the  influence  coefficient  vector  listing  the  displecement  of  the  degrees  of  freedom  due  to 
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a  unit  support  displacement,  and  e(< )  is  the  acceleration  history  applied  at  the  support. 


Dynamic  Response  Analysis  Procedure 

The  orthogonality  properties  identified  in  the  derivation  of  the  Lanczos  algorithm  can  be 
utilized  effectively  in  transforming  the  equations  of  motion  to  a  reduced  form.  Equating  the  two 
expression  for  r;  given  by  Eqs.  16b  and  16e,  and  substituting  Eq.  16a  leads  to 

r,  —  0,+i  q,+i  -  K_l  M  q,  -  q,  a,  -  q,_,  0,  (19) 

Transferring  the  right  side  to  the  left,  the  combined  set  of  the  Lanczos  vectors  obtained  after  m 
steps  can  be  arranged  in  a  single  equation  as  follows: 


> 

•  1 

’ 

M 

K'M 

Qm 

- 

Q. 

0 

■ 

.  , 

»  , 

r.  ®*  (20) 


where  Q„  is  the  n  X  m  matrix  with  columns  qj,  <fe,  -  -  -  q„ ,  e„  is  the  last  column  of  the  I„ , 
identity  matrix  of  dimension  m,  and  T„  is  a  tridiagonal  matrix  made  up  of  the  coefficient  a,,0,: 


<*1  02 
02  ®2  03 

03  .  . 


T, 


(21) 


•  •  0- 

0m  a* 

Now  multiplying  Eq.  20  by  Q^M  and  applying  the  orthonormai  properties  of  the  Lanczos  vec¬ 
tors,  Q-Tmq.  -  Im  leads  to 


Q.rMK-'MQi-T1  (22) 

That  is,  the  projection  of  M  K'1  M  on  the  subspace  with  basis  Q„  is  a  tridiagonal  matrix. 

This  tridiagonal  property  will  be  used  advantageously  in  the  dynamic  response  analysis  pro¬ 
cedure.  For  this  purpose,  the  equations  of  motion  (Eq.  1)  are  multiplied  by  MK~*  with  the 
result: 
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M  K-‘  M  u  +  M  K1  C  u  +  M  u  -  M  K1  f  (23) 

Before  accepting  this  expression  it  should  be  noted  that  the  mass  matrix  M  used  in  dynamic 
response  analyses  frequently  is  singular  (because  mass  coefficients  may  not  be  associated  with  all 
degrees  of  freedom),  and  also  that  significant  information  may  be  lost  when  a  set  of  equations  is 
multiplied  by  a  singular  matrix.  Therefore,  the  use  of  Eq.  23  may  not  be  justified.  However,  it  is 
shown  in  the  Appendix  that  using  Eq.  23  instead  of  Eq.  1  results  in  the  same  solution,  whether  or 
not  M  is  singular. 

Now  expressing  the  response  in  Lanczos  coordinates,  u  =  Q„x(<)  and  completing  the 
transformation,  Eq.  23  becomes 

QlM  K  1  M  Q„x  +  K'1  C  Q„x  +  Qj\f  Qmx  =  Q«M  K  l  f 

Neglecting  damping  and  using  Eq.  22  and  the  orthonormal  property  of  the  Lanczos  vectors;  this 
may  be  written  in  tridiagonal  form: 

T.»+x*l.  (24) 

in  which  g„  =  Q*  M  K'1  f.  If  Rayleigh  damping  is  used: 

C  *  +  atK 

the  tridiagonal  form  is  still  preserved,  as  follows 

Tm  x  +  (o0Tm  +  a,Iw)  x  +  x  =  gm  (25) 

If  the  more  general  form  of  proportional  damping  given  by  the  Caughey  series  [8]  is 
adopted: 

C-M  £  m*[m-1k]‘ 

■* 

then  the  corresponding  damping  term  to  be  used  in  Eq.  24  will  be  of  the  form 

C"  “ 

It  should  be  noted  that  the  damping  only  approximates  the  projection  of  the  Caughey  damping 
on  the  Lanczos  subspace  but  the  discrepancy  is  of  no  significance  in  practical  cases. 

The  tridiagonal  for  of  Eq.  24  or  25  is  particularly  advantageous  when  the  applied  load  is  of 
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the  form  of  Eq.  17  because  then  the  forcing  function  can  be  reduced  as  follows: 

g,  -  Q.  M  K1  f 
■  Q,MK',b«(<) 

■Q«Mq xfot  (I) 

=  Q.  M  Q,  «i  A  c  (<) 

-  A  *!«(<)  (26) 

where  e,  is  the  first  column  of  the  identity  matrix  I, .  Eq.  24  then  takes  the  form 

T.!+  x-  A«i7(0  (27) 

which  shows  that  the  excitation  is  applied  only  in  the  first  of  these  equations  of  motion.  Of 
course,  the  same  loading  would  be  applied  to  the  damped  system  of  Eq.  25. 

The  solution  to  the  equations  of  motion  in  the  form  of  Eqs.  24  or  25  may  be  obtained  either 
by  a  time-stepping  scheme  such  as  the  Newmark  method  [10],  or  the  eigenvectors  of  T„  can  be 
calculated  and  the  response  obtained  by  superposition  of  the  uncoupled  equations.  The  storage 
demands  and  the  number  of  operations  required  for  either  of  these  analytical  procedures  is  of  the 
order  m .  This  is  in  contrast  with  the  cost  of  solving  an  r.  x  n  eigenproblem  in  order  to  employ  a 
standard  mode  superposition  analysis.  Moreover,  because  the  Lanczos  vectors  explicitly  recognize 
the  load  distribution,  fewer  Lanczos  vectors  than  mode  shapes  are  needed  to  obtain  a  desired  level 
of  approximation. 

Required  Number  of  Vectors 

In  order  to  determine  how  many  Lanczos  vectors  may  be  required  to  obtain  the  desired 
accuracy  in  a  dynamic  analysis  it  is  necessary  to  assess  the  errors  resulting  from  a  given  approxi¬ 
mation;  or  more  precisely,  to  determine  the  number  of  vectors  required  with  a  tolerable  error. 
The  participation  factor,  which  indicates  the  component  of  the  applied  load  that  contributes  to 
the  response  of  the  corresponding  Lanczos  vector,  provides  a  measure  of  the  significance  of  that 
vector  in  the  total  response.  Defining  the  j  participation  factor  as 

H,  -  q/*»  (28) 

the  product  Qj  b«h,  gives  a  vector  listing  the  participation  factors  of  all  m  Lanczos  vectors. 


•5/e. 
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So  if  the  elements  of  this  vector  are  evaluated  sequentially  during  construction  of  the  Lancsos 
vectors,  their  successive  values  may  be  judged  as  a  criterion  for  termination  of  the  Lanczos  algo¬ 
rithm. 

The  values  of  the  participation  factors  associated  with  each  Lanczos  vector  may  be  obtained 
at  negligible  cost  during  the  determination  of  the  vector.  Multiplying  Eq.  16b  by  br  yields 

br  r,  =  br  7,  -  brq;  a,  -  br  q;_,  fi,  (29) 

Noting  from  Eq.  16a  that  F;  =  K~lMq; ,  and  recalling  that  the  starting  vector  was  given  by 
$,qt  =“  K  ‘b  (Eq.  18),  the  first  term  on  the  right  side  of  Eq.  29  becomes 
brK~lMq;  —  /J1q1TMq;  -=  0.  Thus,  introducing  Eq.  16c  on  the  left  side,  Eq.  29  can  be  reduced 
to 


0,+i  bT  q,+,  —  -  brq,  a,  -  br  qM  fi, 
or.  making  use  of  the  definition  of  the  participation  factors  (Eq.  28),  this  gives 


ty+i 


{ct,n,  +  0,  >1,-1 ) 

0J+1 


(30) 


Including  this  simple  scalar  calculation  with  the  Lanczos  algorithm  makes  it  possible  to  terminate 
the  calculation  when  the  magnitude  of  »j,+1  drops  bellow  a  specified  tolerance  value. 


Loan  of  orthogonality 

The  simple  Lanczos  algorithm  presented  here,  involving  orthogonalization  with  only  the  two 
preceding  vectors  at  each  step,  is  subject  to  loss  of  orthogonality  with  respect  to  earlier  vectors 
due  to  roundoff  error.  If  such  errors  are  not  corrected  when  they  reach  a  critical  size,  the  vectors 
may  become  linearly  dependent.  Paige  [11]  characterized  the  way  in  which  orthogonality  is  lost 
and  provided  the  theoretical  background  for  avoiding  this  gradual  loss  of  orthogonality  without 
applying  full  reorthogonalization  at  each  step.  From  Paige's  work  a  number  of  procedures  were 
developed  (1,12]  which  take  advantage  of  the  simplicity  of  the  basic  Lanczos  algorithm  while 
maintaining  robust  orthogonality. 

An  algorithm  proposed  by  Simon  (12|  is  adopted  here  as  the  simplest  means  of  retaining  the 
necessary  level  of  orthogonality.  The  vector  w;  —  q;iflMQ;  is  formed  to  indicate  the  loss  of 
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orthogonality  at  may  step  j;  that  is  w,  contain*  tenon  of  the  order  of  roundoff,  aad  Simon's  algo* 
rithm  in  aa  inexpensive  scheme  for  updating  w,  at  each  step  of  the  Laacson  analysis.  His  formu¬ 
lation  for  w;+t  may  be  written 

4,+  i  w,+I  ■  T,  w;  -  *,  w,  -  0,  wM  (31) 

To  start  thin  algorithm  it  is  necessary  to  note  that  w0  —  0  aad  wt  *•  qj/Hqj;  then  the  farther 
values  of  w;  are  obtained  at  the  cost  of  only  5  j  multiplications  per  step.  The  magnitude  of  the 
elements  of  this  vector  demonstrate  the  trend  of  the  orthogonality  condition.  When  any  element 
of  w;  is  excessive,  q,+i  is  orthogonaliied  against  all  the  preceding  vectors.  Then  the  simple 
Lancsos  algorithm  is  continued  until  the  vector  w,  at  some  subsequent  step  again  indicates  the 
need  for  orthogonalization  against  all  vectors. 


Numerical  Results 

As  an  example  we  use  the  finite  element  model  of  a  multistory  building  frame  shown  in  Fig. 
1.  The  model  consist  of  170  two  dimensional  truss  elements  resulting  in  a  total  of  100  degree’s  of 
freedom.  The  material  properties  were  chosen  such  that  the  fundamental  period  is  close  to  1 
second  (Tj  »  1.082  sec.).  A  step-function  load  is  applied  at  the  third  story  of  the  structure  as 
shown  in  the  figure.  This  loading  configuration  ensures  the  participation  of  the  higher  modes  in 
the  response  of  the  structure.  It  is  worth  noting  that  this  model  also  may  be  viewed  as  a  lattice 
analogy  representation  of  a  plane  stress  orthotropic  elasticity  problem,  consequently  the  response 
to  the  step  function  load  will  contain  the  effects  of  a  stress  wave  rediating  from  the  load  point. 
Also  it  should  be  noted  that  the  structure  was  assumed  to  be  undamped  in  order  that  the  higher 
mode  contributions  to  the  response  would  not  be  suppressed. 

The  Lanczoe  method  was  applied  to  reduce  the  equations  of  dynamic  equilibrium  to  the  tri- 
diagonal  form.  The  Newmark  step-by-step  method  then  was  used  to  evaluate  the  response  of  the 
reduced  problem.  Newmark  parameters,  0  and  7  were  chosen  to  be  0  ■■  •£•  and  7  m  1  to  avoid 
any  numerical  damping  for  consistency  with  the  physical  assumption  and  the  time  step  was  taken 
to  be  ■■  0.001  sec.  This  requires  performing  just  over  1000  time  steps  per  fundamental  period 


i 


13- 


of  tke  response.  Tke  analysis  was  carried  oet  for  various  Bombers  of  Lancxos  vectors,  la  figure  2 
the  displacement  response  at  tke  top  story,  a,  Is  compared  with  that  obtained  by  applying  the 
Newmark  method  to  the  on  red  need  problem  using  the  same  parameters.  From  figure  2  it  is  clear 
that  the  response  of  the  structure  can  be  obtained  accurately  with  as  few  as  IS  Lancxos  vectors. 
In  figure  3  we  show  the  way  tke  displacement  «  at  time  t  ■  3.0  sec.  stabilises  as  the  number  of 
Lancxoa  vectors  is  increased.  This  convergence  behaviour  is  also  typical  of  the  finite  element 
method  as  tke  number  of  elements  are  increased.  This  similarity  is  dae  to  the  fact  that,  like  tke 
finite  element  method,  the  Lancxoa  method  also  minimises  the  potential  energy.  From  this  plot  it 
can  be  seen  that  with  only  6  Lancxoa  vectors  the  displacement  u  is  to  within  5 %  and  with  20  vec¬ 
tors  to  within  2%  of  tke  solution  obtained  by  the  direct  method.  Figure  4  shown  the  variation  of 
the  Lancxoa  participation  factors,  iy, ,  with  tke  number  of  Lancxoa  vectors.  It  is  apparent  that  aa 
rjj  approaches  xero  the  displacement  response  stabalixes.  Therefore  as  soon  as  if,  falls  bellow  a 
specified  tolarence  the  Lancxos  algorithm  can  be  terminated. 

The  orthogonality  state  among  the  Lancxos  vectors  was  monitored  using  the  algorithm  of 
Eq.  31.  In  this  example,  a  reorthogonali  ration  step  was  carried  out,  on  average,  every  5  iterations 
with  the  first  occuring  at  iteration  <J. 

Concluding  Remarks 

This  paper  points  out  that  the  Lancxos  algorithm  can  be  related  directly  to  the  inverse 
iteration  method  of  eigenvector  evaluation.  The  essential  concept  of  the  Lancxos  procedure  is 
that  the  successive  vectors  obtained  during  inverse  iteration  become  the  coordinates  used  in  the 
dynamic  response  analysis  when  they  are  successively  made  orthogonal  to  the  previously  derived 
vectors.  It  is  shown  here  that  the  orthogonality  need  be  enforced  only  with  respect  to  the  two 
preceding  vectors;  the  way  the  Gram-Schmidt  orthogonalixation  is  applied  automatically  easurea 
orthogonality  with  all  earlier  vectors,  and  this  fact  permits  the  equations  of  motion  to  be  restated 
in  tridiagonal  form.  Thus  the  dynamic  respoase  may  be  obtained  by  step-by-step  integration  of 
the  simplified  equations. 


Aa  ittnntru  mljrtii  procedure,  of  com—,  would  be  to  derive  tke  ribntki  mode  ship— 
from  the  Laactos  coordinates,  — d  (lea  to  calculate  the  teapooae  ia  ten—  of  the—  ntoaphd 
nodal  coofdmatea.  However,  boeaa—  the  dyaaaric  loadtag  geaeraly  ia  defied  —  a  aomerical 
aeq—ce,  the  reapo—e  usually  mast  be  ohtaiaed  by  step-by-etep  method  evea  for  a  mode  superpo¬ 
sition  notation.  Typical  time' stepping  methoda  a— h  —  tho—  of  the  Newmark  family,  requires  a 
total  of  8m  +  (4  to  -  2)m  operations  per  time  step  to  evaluate  the  eolation  to  aa  —damped  sys¬ 
tem  of  m  eqaationa  with  its  matrices  haviag  a  half  band  width  w.  When  the  system  is  tridiago- 
aal  (v  «■  2),  aa  ia  Eq.  27,  the  cost  per  time  step  becomes  14m  operations,  aid  when  the  system  is 
u coupled  this  cost  oily  reduces  to  10m.  It  appears  that  this  additional  redaction  is  not  substan¬ 
tial  enough  to  justify  the  overhead  cost  of  3m*  +  9m*  operations  required  for  reducing  Eq.  27  to 
aa  —coupled  system  of  equations. 
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Appendix 

The  question  that  is  addressed  here  is  whether  the  solution  to  Gq.  23  is  equivalent  to  the 
solution  of  the  original  equation  (Eq.  1)  when  M  is  singular.  For  simplicity  the  undamped  case  is 
considered,  and  with  no  loss  in  generality,  the  mass  and  stiffness  matrices  can  be  partitioned  as 
follows: 


FAX  0]  fiCu  K12 

0  oj  “d  K  “  (k2i  Kas 

where  Ku  =  Ka-  The  explicit  forms  of  K  and  M  is  used  to  reduce  Eq.  1  to 


M]  U]  +  Kjj  U]  +  Ku  u2  •  ffi 
+  Ka  #i  *■'  *a 


(32a) 

(32b) 


where  u  *■  j”‘ j  and  f  *■  JJ'j.  Multiplying  Eq.  32b  by  K12Ka  and  subtracting  the  result  from 
Eq.  32a,  eliminates  u2  to  get 


Mt  Ui  +  ( Ku  -  K„  K&  Ka  )  Hi  “  fi  -  Ku  Ka  fs  (33) 

The  next  step  is  to  obtain  an  expresion  for  K'1  by  performing  the  block  Lr  DL  fac  tori  ra¬ 
tion  of  K  (not  the  usual  LDLr)  and  inverting  the  result  to  get 


r-l . 


I  0 

,-K^Ka  I 

(Ku-K^K ar*  0  1  |l-Kl2Ki* 


0  tU 


(34) 


When  putting  this  relation  and  the  singular  expresion  for  M  in  Eq.  13,  the  first  term  on  the  left 
hand  side  becomes 


M  Kl  M  u 


M,  (Kn-Ki8K5j  Ka M  i  St 


(35) 


and  similarly  for  the  right  hand  side  of  Eq.  13  yields 


MK‘f- 
and  the  first  of  these  equations  becomes 


M,(K„  -  KuKa  Kn )  (f,  -  K12Ka  f2) 
0 


(38) 


f 


M,  (Kn-KuKaKjjf  MjUi  +  Mj«t  -  Ml[Kn~KuK^Kn)1  (f.-K^K^f,) 
Mult  Ip  lying  this  equation  by  (Ku-Ku  KgKa)  Mf1  reduces  to  Eq.  33. 
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Young's  Modulus  »  80.0  MN/ram* 
Croes-sectionsl  Ann  »  0.01  m1 

Mass  Density  —  l8.0X10*kg/m1 


Figure  1.  Details  of  the  Truss  Example. 
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